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We study localized nonlinear excitations of a dilute Bose-Einstein condensate (BEC) with spin- 
orbit coupling in a deep optical lattice (OL). We use Wannier functions to derive a tight-binding 
model that includes the spin-orbit coupling (SOC) at the discrete level in the form of a generalized 
discrete nonlinear Schodinger equation. Spectral properties are investigated and the existence and 
stability of discrete solitons and breathers with different symmetry properties with respect to the 
OL is demonstrated. We show that the symmetry of the modes can be changed from on-site to 
inter-site and to asymmetric modes simply by changing the interspecies interaction. Asymmetric 
modes appear to be novel modes intrinsic of the SOC. 

PACS numbers: 


Introduction. Presently there is a growing interest in 
the study of Bose-Einstein condensates (BEC) in the 
presence of non abelian gauge fields that mimic magnetic 
interactions. In particular, spin-orbit couplings (SOC) 
of Rashba and Dresselhaus types and their combina¬ 
tions, have been recently realized in binary BEC mixtures 
in thcpresence of trapping potentials V{x) of different 
types Hi (see the review articles @, 0 )- As well known, the 
SOC represents a major source of intra-atomic magnetic 
interaction. In solid state physics it plays an important 
role in the magnetism of solids, well described in terms of 
individual ions, as it is for earth rare insulators. In gen¬ 
eral, however, in solids the SOC is a rather weak source 
of magnetic interaction (largely superseded by the elec¬ 
trostatic effects, which is impossible to manage/enhance. 

The situation is quite different in Bose-Einstein con¬ 
densates where a variety of forms of synthetic spin-orbit 
couplings can be easily generated by external laser fields 
and the strength can be easily controlled. BEC with 
SOC, indeed, now represent ideal systems to explore in¬ 
teresting phenomena that are difficult to achieve in solid 
state, such as new quantum phases with unusual mag¬ 
netic properties, existence of stripe modes 3 , fractional 
topological insulators@-0|, Majorana fermions, etc. 

The interplay of the spin-orbit coupling with the inter¬ 
atomic interactions (nonlinearity) and the periodicity of 
the OL also leads to the existence of SOC solitons[ 8 |, Iff. 
Gap solitons in BEC with periodic Zeeman field were re¬ 
cently discussed in 10] (see also the review of BEC in 
OL with SOC in ll[). In particular, different settings 
for the spatial periodicity have been considered: peri¬ 
odicity in each separate component HU. [Til ], periodicity 
in the Raman coupling fl4|. periodicity in the Zeeman 
field E Dispersion relations of one-dimensional 

BEC with SOC in OL were experimentally investigated in 
fist and existence of flat bands and superfluidity in BEC 
with SOC demonstrated in Array of vortex lattices 


in shallow optical lattices [l3| and vortices in 2D opti¬ 
cal lattices jl9[ were also investigated. The extension to 
quasi-periodic OL and the Anderson localization of BEC 
with effective SOC have been considered in {20|. These 
studies refer mainly to the continuous case, and the ef¬ 
fects of SOC on lattice models (BEC arrays) is practically 
uninvestigated. 

The aim of the present letter is to consider BEC mix¬ 
tures with SOC in the presence of deep optical lattices. 
We use Wannier functions to derive a tight-binding model 
that includes the SOC at the discrete level in the form 
of a generalized discrete nonlinear Schodinger equation 
(SO-DNLS) with typical double minima linear dispersion 
relations, supporting nonlinear localized two component 
excitations with chemical potentials inside the forbidden 
zone of the band structure (gap-solitons), in the non¬ 
linear case. The existence and stability of gap-solitons 
(GS) has been investigated as a function of parameters 
such as the strength of the interatomic interactions and 
strength of the SO coupling. The band structure and the 
location of the modes in the band structure is shown to 
be symmetric with respect to the change of the signs of 
both inter- and intra- species interactions. We find that 
for attractive (resp. repulsive) nonlinear interactions, GS 
can exist in the semi-infinite gap below (resp. above) the 
bottom (resp. top) band, as well as in the gap between 
the two bands. The time evolution of attractive (resp. 
repulsive) fundamental GSs, e.g. the one with the lower 
(resp. higher) chemical potential, are very stable while 
the ones located in the intraband gap are typically un¬ 
stable. Quite interestingly, by increasing the strength of 
the intra-species interaction 7 , (for simplicity assumed 
equal for both species) away from the linear limit, for a 
fixed and equal sign of inter-species nonlinearity, we find 
three distinctive regions in which GS undergoes sponta¬ 
neously symmetry breaking. More precisely, in the range 
0 < |y| < | 7 i| the GS are found to be asymmetric with re- 
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spect to the lattice points, in the interval |7i| < I7I < I721 
they display inter-site symmetry (eg symmetry with re¬ 
spect to the middle point between two consecutive lattice 
sites) and above I72 1 display the on-site symmetry, this 
behavior being observed both for attractive and repulsive 
case. Asymmetric modes appear as novel modes induced 
by the SOC. 

A possible physical interpretation of this phenomenon 
is suggested which could be valid also for the correspond¬ 
ing continuous case. 

Tight-Binding Model. We consider BEC with equal 
contributions of Rashba and Dresselhaus SOC that can 
be described in the mean field approximation by the fol¬ 
lowing coupled Gross-Pitaevskii equations [ 1 ]: 



(-^2 + V ( x ))^o - + 

(9j\^j \ 2 + j = 1,2, (l) 


with the linear coefficients a, O = ffi = — O2 arising 
from the spin orbit interaction while the nonlinear ones, 
g and gi,i = 1,2, are related to inter-species and intra¬ 
species scattering lengths, respectively. In the following 
we consider as trapping potential an optical lattice (OL), 
e.g. a periodic potential of the form V(x) = Vqcos{2x) 
(equal for the two components), and concentrate on the 
case of large amplitudes Vo >> 1 (deep optical lat¬ 
tice) for which it is possible to develop the tight bind¬ 
ing approximation[ 21 J. To this regard we expand the 
two component fields in terms of the Wannier functions 
w{x — n) localized around lattice sites n of the underlying 
uncoupled linear periodic eigenvalue problem, e.g. with 
a = Cl = 0 ,g = gi = 0 , in Eq. ([[]), 


i 

2 "b ^(*^) Tm,k = £m{k)<Prn,k- (2) 


Here g} m ,k and £ m {k) denote, as usual, Bloch (Floquet) 
functions and energy band, with m the band index and k 
the Bloch wavenumber taken in the first Brillouin zone. 
For a deep optical lattice it is convenient to use the Wan¬ 
nier basis and expand the fields as 


4 >i 


5 ~ 2 u n (t)w m {x-n ), t/>2 = '^ 2 ,Vn(t)w m {x-n), ( 3 ) 

n,m n,m 


where w m (x — n) denote Wannier functions associated to 
the periodic problem © , and u n , v n are time dependent 
expansion coefficients to be fixed in such a manner that 
Eq. (HD are satisfied. Notice that in the expansions 0 
appear the same Wannier functions since the underlying 
periodic problem is the same for both components. Also, 
we remark that Wannier functions are orthonormal with 
respect to both the band index and the lattice site n 
around which they are centered: 


J Wj(x — n)*Wji(x — n')dx = Sjj'dn^i ( 4 ) 


In the following we restrict only to one band (the ground 
state band) and therefore we drop out the band index 
m. By substituting 0 in Eq. © and projecting the 
resulting equations along the w(x — n ) function, one ar¬ 
rives at the following coupled SO- DNLS system for the 
coefficients 


= ~ r ( Mrl + 1 + U n-l) + - V n -!) + Llu n + 


(ll\ u n\ 2 + r y\v n \ 2 )u n , 
Vn+1 + V n -1 
{l\u n \ 2 +l2\v n \ 2 )Vm 


= - r(v n+ l + V n - 1) + i — (u n+ l - U n -1 ) - + 


( 5 ) 


with 


/' d 2 

r = Tn t n +1 = / w(x - n)*—w(x - (n + 1), 

7 i= 9 i J \w(x - n)\ 4 dx, 7 = g j\w{x - n^dx, (6) 

f d 

a = ain , n + 1) = 2a / w(x — n)*—w(x — (n + 1)). 
J ox 

In the derivation of Eq. © the following relations among 
coefficients have been used: 

er(n, n) = 0, er(n, n — 1) = — <r(n — 1, n) = —a (n, n + 1), 


r(n, n + 1) = T(n, n — 1). 


Moreover, due to the strong localization of the Wannier 
functions around the lattice sites, the sums on n have 
been restricted to onsite and to next neighbor sites only, 
for non diagonal and diagonal terms, respectively. 

Linear Case. Let us first consider the case with all 
nonlinear coefficients detuned to zero, 7 = 7* = 0, 


. du n 

i—r~ 

dt 

. dVn 
dt 


r(tin-(-l —t Un — l) T i 2 (^71+1 Un—l) T ^lUniC) 
r(Un+l ~h Vn — l) T i-^ (Un+l Un — 1) HUn-(8) 


corresponding to the case of a linear chain with spin- 
orbit coupling. From physical point of view this corre¬ 
spond to a one-dimensional array of BEC in deep optical 
lattice in the presence of SO interaction. To find the dis¬ 
persion relation of this linear chain we consider solutions 
of the form 


_ A^i(kna—ujt) ^ ^i(kna-Lat) 


( 9 ) 


with wavenumber k n — 2 n/Ln varying in the first Bril¬ 
louin zone [—7r/a,7r/a], with —N /2 <n< N/ 2 , L = Na 
the length of the chain, a the lattice constant fixed below 
without loss of generality to a = 1 . Substituting © in 
Eq. © one obtains a homogeneous system of equations 
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FIG. 1: Chemical potential vs k in the reciprocal space for the 
DNLS chain with SO coupling for the linear case 71 = 72 = 
7 = 0 and parameter values T = 0.3, fi = 1.352, a = 1.5. The 
red and blue dots represent chemical potentials of localized 
modes for the nonlinear cases: 71 = 72 = —0.65, 712 = —1.8 
(attractive case), and 7 ! = 72 = 0.65, 712 = 1.8 (repulsive 
case), respectively. Modes are plotted on the k = 0 line for 
graphical convenience. 


for the coefficients A, B , whose compatibility conditions 
directly leads to the dispersion relation 

uj{k)± = —2r cos (k) ± \JQ 2 + a 2 sin 2 (fc). (10) 


It is worth to note that the upper and lower branches 
of the dispersion curve are related as follows 


W±(fc) = —Wzp(fc + 7r). (11) 


Also notice the presence of two degenerated minima at 
±fc_ on the ground state branch and two degenerated 
maxima at ±fc + on the upper branch, with k± given by 


k± = ±arcos 


2T I a 2 + n 2 I 

a V <t 2 + 4r 2 . ' 


( 12 ) 


The dispersion curve in the first Brillouin zone is depicted 
for typical parameter values in Fig. |T| One can readily 
show that the amplitudes of the two components must 
be related by 


(f )± = T y/n 2 + <T* sin \k)). (13) 


By combining the k + and fc_ modes one readily obtain 
stationary stripe solutions for w+ (top panels of Fig. [2J 
and w_ (bottom panels of Fig[5j). 

Nonlinear Case. In the presence of the nonlinearity 
the existence of localized modes with chemical potentials 
in the forbidden zones of the bandgap structure (discrete 
gap-solitons) become possible. 

In particular, in analogy with gap solitons or intrin¬ 
sic localized modes of continuous and discrete Gross- 
Pitaevskii equations in absence of spin orbit coupling, 
one expects that for attractive (resp. repulsive) inter 
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FIG. 2: Panel (a). Existence curves of the localized modes 
of the DNLS with SO coupling as a function of the equally 
attractive intraspecies interaction 7 = 71 = 72 for parameters 
T = 0.3,12 = 1.352, a = 1.5 and fixed intraspecies interaction 
712 = —1.8. The inset display details of the bottom curves, 
enlarged along the vertical axis. Panels (b), c and (d) show 
Un (blue lines) and v n (red lines) profiles of gap solitons at 
7 = —0.65 (panel (b)), 7 = —0.52 (panel (c)) and 7 = —0.3 
(panel (d)), corresponding to red points (bottom profiles), 
blue points (middle profiles) and black points (top profiles) 
shown in panel (a) (in all cases v n is purely real while u n 
is purely imaginary). Middle and top profiles in all panels 
are shifted upward by 1.5 and 2.5, respectively, to avoid over¬ 
lappings. Vertical dotted lines separate regions of different 
symmetry type. All plotted quantities are in dimensionless 
units. 


and intra-species interatomic interactions, gap-solitons 
to form with chemical potentials just below local min¬ 
ima (maxima) of the dispersion relation. This is indeed 
what one finds from a numerical self-consistent diagonal- 
ization of the stationary eigenvalue problem associated 
to Eq. (HI) (see Fig. JT]) for typical examples). 

In the panel a) of Fig. [2] are shown existence curves 
of gap-solitons of the DNLS with SO coupling for the 
case 712 < 0 and for equally attractive intra-species in¬ 
teractions 71 = 72 = 7 < 0 . The lower two branches 
correspond to the ground localized modes in the lower 
semi-infinite gap while the top curve refers to a mode in¬ 
side the interband gap. In remaining panels of the figure 
are shown the imaginary and real parts of the component 
profiles u n (blue lines) and v n (red lines) of the gap soli¬ 
tons for the different values of 7 corresponding to black, 
blue and red points depicted in panel (a). 

It is worth to note that there is a phase difference of 
e J 7 r /2 between the u n and v n components (v n being real 
and u n purely imaginary) as well as different symmetry 
properties with respect to the lattice sites for the two 
components. Also notice that the symmetry properties 
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FIG. 3 : Time evolution of the localized discrete modes shown 
in panels b of Fig. [ 2 ] to the parameter value 7 = — 0 . 65 .Top 
and bottom panels refer to the ground state and to the first ex¬ 
cited state in the lower semi-infinite gap, respectively. Other 
parameters are fixed as in Fig. [ 2 ] 
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FIG. 4 : Time evolution of the localized discrete modes shown 
in the panel c of Fig. [ 2 ] for the case 7 = — 0 . 52 . Top and 
bottom panels refer to the ground state and to the first ex¬ 
cited state in the lower semi-infinite gap, respectively. Other 
parameters are fixed as in Fig. [ 2 ] 


FIG. 5 : Time evolution of the localized discrete modes shown 
in the panel d of Fig. [ 2 ]for the case 7 = — 0 . 3 . Top and bottom 
panels refer to the ground state and to the first excited state 
in the lower semi-infinite gap, respectively. Other parameters 
are fixed as in Fig. [ 2 ] 

of the modes change as the intra-species interactions are 
varied. In general, the following situation is observed. 
For a fixed attractive inter-species nonlinearity, we find 
three distinctive regions in which GS undergoes sponta¬ 
neously symmetry breaking as the strength of the attrac¬ 
tive intra-species interactions (assumed equal e.g. 7 i = 7 , 
for simplicity) away from the 7 = 0 limit. More precisely, 
in the range —0.35 < 7 < 0 the GS are found to be asym¬ 
metric with respect to the lattice points, in the interval 
—0.6 < |y| < —0.35 they display a symmetry with re¬ 
spect to the middle point between two consecutive lattice 
sites (inter-site symmetry) and for 7 < 0.6 they display 
the on-site symmetry. Notice that the borders of these 
regions correspond to the appearance of small kinks in 
the existence curves /r versus 7 and have been evidenced 
by dashed vertical lines (see inset in panel a) of and Fig. 
El- Also notice that the kinks at the left border of the 
middle region separating the modes with on-site symme¬ 
try (see panel (b)) from the ones with intra-site symme¬ 
try (see panel (c)), both of symmetric or anti-symmetric 
type, are sharper than the ones on the right border where 
modes becoming asymmetric (see panel (c)). This also 
correlates with the fact that in general one would expect 
that a lost of symmetry being a smoother process with 
respect to a spontaneously change of a symmetry. 

The dynamical properties of the nonlinear modes de¬ 
picted in Fig. Elc)-(d)) have been investigated by direct 
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FIG. 6 : Time evolution of the the on-site discrete gap soliton 
inside the inter-band gap shown in the panel b of Fig. [2] (top 
profiles) for the case 7 = —0.65. Other parameters are fixed 
as in Fig. [2] 


numerical integrations of the DNLS-SO, by taking them 
as initial conditions with a small noise component im¬ 
posed to check their stability under time evolution. This 
is shown in Figs. [3][5]for the two lowest nodes in the semi¬ 
infinite gap. As one can see, for the chosen parameters 
the on-site symmetric modes are both very stable (see 
Fig. 0 while for inter-site symmetric modes stability is 
achieved only for the ground state, the first excited state 
being metastable. Notice, in the last case, that the mode 
decays into an on-site- symmetric mode plus background 
radiation, similarly to the inter-site-on-site transition ob¬ 
served for discrete breathers of coupled DNLS in absence 
of SO coupling. Quite surprisingly, asymmetric modes 
of the region 7 < —0.3 also appear to be stable (or long 
lived) under time evolution, as one can see form Fig. [5] 
The existence of such modes appears to be typical of the 
SO coupling, since they would not be possible just in 
ordinary coupled DNLS. The stability properties of lo¬ 
calized modes in the semi-infinite gap, however, appear 
to be critical. In particular, for the parameter ranges we 
have explored and independently from their symmetry 
types, none of them was found to be stable (see Fig[6]for 
a typical example). 

Similar results are found for the case of all repul¬ 
sive interactions. By reversing signs of all the interac¬ 
tions chemical potentials of the discrete gap solitons also 
change their signs (see Fig. [[]). Shape and symmetries 
of the modes, however, are different in the two cases as 
one can see by comparing Figs. [7] with [2j Also notice 
the symmetry of corresponding existence curves in the 
p — 7 plane, with appearance of kinks at the change of 
symmetry points. The stability properties under time 
evolution, look also similar to the all attractive case. In 
particular, modes inside the inter-band gap are also found 
to be unstable, while the ones with highest chemical po¬ 
tentials inside the upper semi-infinite gap are typically 
stable (not shown for brevity). 

In closing this Letter we provide parameter estimates 
for possible experimental implementations of the above 
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FIG. 7: Panel (a). Existence curves of the localized modes 
of the DNLS with SO coupling as a function of the equally 
repulsive intraspecies interaction 7 = 71 = 72 for parameters 
7 = 0.3, Q. = 1.352, <7 = 1.5 and fixed intraspecies interaction 
71.2 = 1.8. The inset display details of the top curves, enlarged 
along the vertical axis. Panels (b), (c) and (d) show u n (blue 
lines) and v„ (red lines) profiles at 7 = 0.65 (panel (b)), 7 = 
0.52 (panel (c)) and 7 = 0.3 (panel (d)), corresponding to 
black points (bottom profiles), blue points (middle profiles) 
and red points (top profiles) in panel (a), respectively (in all 
cases v n is purely real while u n is purely imaginary). Middle 
and top profiles in the b,c,d, panels are shifted upward by 1.5 
and 2.5, respectively, to avoid overlapping. Vertical dotted 
lines separate regions of different symmetry type. All plotted 
quantities are in dimensionless units. 


results. For this one can adopt an experimental setting 
similar to the one in jl|. In particular, the SOC can be 
realized in 87 Rb using two counter-propagating Raman 
lasers with Xr = 804nm, spin-orbit coupling a = hk^/m , 
and hfl = 2Er (here Er = h 2 kj t /2m is the recoil energy). 
The optical lattice can be generated by means of two 
additional beams with Xl = 1540nm[jji]. Typical values 
of the OL potential strength Vo suitable for the tight- 
binding case we considered can be Vo/E r >10 [Illli. In 
the case of all attractive interspecies interaction nonlinear 
modes below the bands should appear for wide range 
of number of atoms, being typically stable. Change of 
the interspecies scattering length to observe the change 
of symmetry of the discussed modes can be achieve by 
means of the Feshbach resonance technique. 


In conclusion we have derived a tight-binding model 
for BEC with SOC in deep optical lattices and demon¬ 
strated for this model the existence and stability of differ¬ 
ent types of discrete solitons. We showed that nonlinear 
modes can change symmetry from on-site symmetric to 
inter-site symmetric and to fully asymmetric as the inter¬ 
species interaction is varied. Asymmetric modes appear 
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to be intrinsic novel excitations of the BEC with SOC. 
The critical values at which the symmetry changes are 
found to correspond to kink-like profiles in the chemical 
potential existence curves that suggest the occurrence of 
phase transitions. The possibility to observe these phe¬ 
nomena in real experiments was suggested. 
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